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Abstract Understanding the reasons of the cyclic vari- 
ation of the total solar irradiance is one of the most chal- 
lenging targets of modern astrophysics. These studies 
prove to be essential also for a more climatologic is- 
sue, associated to the global warming. Any attempt 
to determine the solar components of this phenomenon 
must include the effects of the magnetic field, whose 
strength and shape in the solar interior are far from be- 
ing completely known. Modelling the presence and the 
effects of a magnetic field requires a 2D approach, since 
the assumption of radial symmetry is too limiting for 
this topic. We present the structure of a 2D evolution 
code that was purposely designed for this scope; rota- 
tion, magnetic field and turbulence can be taken into 
account. Some preliminary results are presented and 
commented. 

Keywords Sun: evolution — Sun: interior 



1 Introduction 

It is well known that the Total Solar Irradiance (TSI) 
varies with time, with variations both on a short 
time scale, of the order of days to months, and on a 
long time scale, which fo llows the 11-years solar cycle 
(|Willson fc Hudson! 1 1 9 9 lh . While the short-term vari- 
ations are associated to the solar rotation, as a con- 
sequence of surface spots and faculae, which determine 
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the variation of the energy flux in the Earth's direction, 
the reasons for the long-term variation are more de- 
bated, as confirmed by the extensive literature in that 
regard, ranging from a more "surface" interpretation 
( Krivova et al.l 120031 ) to papers stressing the role of a 
possible variation of the internal structure of the Sun 
( Li fc SofialEioih . 

In this latter case we would expect that at least 
the whole convective envelope is involved, thus opening 
the possibility of variability on even longer time scales 
which could drive long-term variations of the Earth's 
climate on the Earth. This confirms the importance of 
understanding the role played by the internal variation 
of the solar structure. 

Among all the possible causes that may potentially 
alter periodically the solar structure, the magnetic field 
is generally regarded as the main contributor. Some 
preliminary investigations, limited to a ID treatment, 
and thus limited to purely radial symmetry configura- 
tions, showed that the observed variation of the TSI 
can be reproduced by assuming the presence of a large- 
scale magnetic field, which varie s in phase with the solar 
cycle (see e.g. iLi fc Sofia! (|200lh ). 

The ID approach unfortunately proves to be of ex- 
tremely limited relevance for such an important inves- 
tigation, because any realistic model of the Sun, ac- 
counting for rotation and magnetic field, is expected 
to develop asymmetries at least of the order of magni- 
tude of the effects that we want to investigate. Also, it 
is practically impossible to find any morphology of the 
magnetic field that is radially symmetric, whatever its 
toroidal or poloidal components are. 

These arguments stimulated our group to take a step 
further, with the ext ension of the classic YR EC code 
for stellar evolution ( Demarque et al. 2008) to a 2D 



structure, accounting for variations of the main physi- 
cal and chemical quantities not only with the distance 
from the center, but also with the latitude. In this 
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paper we describe the theoretical framework that in- 
troduces bidimensionality, and that accounts for the ef- 
fects of magnetic fields, turbulence and rotation; we 
eventually provide a new set of differential equations 
that must be solved to determine the solar structure, 
and the changes expected when the above mentioned 
phenomena are considered. 

The paper is organized as follows: Sections 2 and 3 
review the modifications to the standard stellar struc- 
ture equations needed to involve a magnetic field and 
turbulence; the structure of the 2D code, with the rele- 
vant equations, is presented in Section 4. Some prelim- 
inary results are given in Section 5. 



2 The inclusion of magnetic field in stellar 
modeling 

The effects due to the presence of an internal magnetic 
field are mode led following the method described in 
Lvdon fc Sofia ( 19951 ). The magnetic field is described 
by means of an energy density \ m — B 2 /8np. We list 
hereafter all the modifications demanded to properly 
account for the presence of this perturbation. 



1. The total pressure is 



Pr = P 



gas 



B 2 

Pad + ^~ 
87T 



where the magnetic contribution was added to the 
standard gas and radiation components. 

2. The equation of state (EOS) is consequently changed,| 
as density depends not only on pressure and temper- 
ature, but also on the strength of the magnetic field: 

d In p = ad In Pt — Sd In T — v m d In Xm 

where a — din p/dhiPx, 8 = dlnp/dhiT, 
v m = dhip/d\nxm 

3. The momentum equation must include a magnetic 
term expressing the Lorentz force, that is 

i-(VxB)xB = -V.(g) + ^.V)S 

Regarding the anisotropic pressure term 

H = V • (B ■ B)/(4np), the latter equation becomes: 



1 

-in 



(V xB)x B = -V(P m ) + pn 



4. The variation of the total entropy St, in agreement 
with the first law of Thermodynamics, must include 
not only the variation of the internal energy and the 



work done by the gas to expand, but also the varia- 
tion of the magnetic energy density. Expressed dif- 
ferentially, the law giving dSx becomes 

TdS T = dU + PdV + d Xm 

5. The adiabatic gradient V a d, defined as the loga- 
rithmic derivative of temperature with respect to 
pressure along a transformation at constant entropy, 
must change, since the entropy itself contains a new 
magnetic term. The complete expression is 



ad 



pc P 



where V m = d In Xm/d In P. 
6. The change in the EOS determines an alteration of 
the criterion for radiative stability, because the lat- 
ter is based on the evaluation of buoyancy acting 
on a convective eddy. The new formulation of the 
Schwarzschild criterion is 



V rad < V ad 1 - — V t 
a 



We note that the presence of a magnetic field with a 
strong gradient in the proximity of the formal border 
of convection may potentially alter the location of 
the convection/radiation boundary. 
The expression for the convective flux is changed in 
the presence of a magnetic field. Within the con- 
text of the Mix ing Length The ory (MLT) modelling 



of convection (jVitense 
Fc 



1953), 



the convective flux 



pvconvdQ is evaluated on the basis of dQ, i.e. 
the excess heat of the convective eddy with respect 
to the environment before dissolving: according to 
expression for the change of the total entropy given 
in point 4 above, inclusion of possible differences in 
the magnetic energy density is required. A further 
modification to the calculation of the convectice flux 
is needed in the evaluation of the convective velocity, 
made on the basis of the work done by buoyancy: in 
this case the excess density is not only related to the 
excess temperature, but also on the different Xm be- 
tween the convective element and the environment 
(see point 2 above). 



The above changes were used bv lLi fc Sofial (|200ll ) 
to perform a qualitative study of the effects of the 
presence of a large scale magnetic field, allowed 
to change periodically with the solar cycle, on the 
global changes of the solar structure. They could 
reproduce the observed 0.1% variation of the so- 
lar luminosity with various magnetic fields, differ- 
ing in intensity and in the region where they peak. 
The main results of this first investigation was that 
deeper magnetic fields required higher intensity, and 
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were associated to stronger variation of the total ra- 
dius of the Sun, than shallower ones. 



3 The role of turbulence 

From the numerical point of view, the inclusion of the 
effects of turbulence is qualitatively similar to consid- 
ering the presence of a magnetic field. The key quan- 
tities in this context are the turbulent kinetic energy 
Xt — v 2 /2, and the degree of anisotropy, commonly ex- 
pressed via the variable j t = l + 2(v r /v) 2 . The presence 
of turbulence also adds a new term to the pressure, i.e. 

Pt = p{lt - i)x*- 

The modifications necessary to include turbulence 
are the same as those given in the previous section, 
provided that Xm and P m are replaced, respectively, by 
Xt and Pt, and that the variation of density in the EOS 
is expressed in terms of dxt and dr/t ■ The magnetic en- 
ergy gradient V m is replaced by V* = dlnxt/dlnP in 
the evaluation of the adiabatic gradient and the super- 
adiabaticity, and in the Schwarzschild criterion. 

The momentum equation is changed accordingly, the 
two additional terms associated to turbulence being an 
isotropic pressure component and an anisotropic term, 
namely 



-V(P t )+pT 



where 



T = V • [p(v r v r - v g ve)e e ee + p{v r v r - v^v^e^e^)/ p 

In the previous expression we have indicated the three 
components of the turbulent velocity with v r , Vg, and 
i>0, while egeg and e^e^ are two of the three components 
of the unit tensor I. 



the same order of magnitude as the perturbations that 
we want to investigate. In addition, the satellite PI- 
CARD, whose launch is schedule for the end of 2009, 
is expected to detect other 2D observational features, 
which can also be used to discriminate among the var- 
ious models. These are the reasons behind our choice 
to develop a full 2D formulation. 

Our approach is based on the equipotential surfaces 
of the gravitational potential, $. The mass coordinate 
is the mass contained within the surfaces $ = const, 
which, unlike the ID case, may not be spherically sym- 
metric. The colatitude 9 is the other independent vari- 
able. The system is assumed to be azimuthally sym- 
metric. Each equipotential surface is characterized by 
the expression r — r($,f?). The density is generally 
not constant on the equipotentials; yet, it is possible to 
introduce an integrated density, defined as 

1 f* 

Pm = J pr 2 sin 9d9 

The basic equations used are the canonical rela- 
tions expressing the equilibrium of the stellar structure, 
namely the conservation of mass, momentum and en- 
ergy, and the Poisson equation. Unlike the ID case, the 
gradient of the gravitational potential $ has a trans- 
verse component Q, and its radial part deviates by 
the unknown quantity U from the standard expression, 
Gm/r 2 : 



(Gm 



u,g 



On the basis of the numerical simulations bv lRobinson 
( 20011 ). which give the valu es of xt and <y t i n the supera- 
diabatic layer of the Sun, iLi et al.l (|2002l ) studied the 
simultaneous effects of magnetic field and turbulence, 
and they presented various models that could repro- 
duce the cyclic variation of the solar irradiance, and 
a nice fit of the observed oscillation frequencies, par- 
ticularly in th e 2000-4000 ^IIz regime (see Fig. 12 in 
Li et al.1 (<2002h b 



The horizontal component of the energy flux, Fg , must 
also be considered, although the luminosity is still re- 
lated to the radial flux. 

In agreement with the discussion of the previous 
c tS£]cJjons, the momentum equation needs modifications 
with respect to the standard case to account for the 
presence of magnetic field and turbulence. Rotation 
generates a centrifugal acceleration that must be in- 
cluded: 1Z = Q 2 r sin 9 (sin 9, cos 9) 

Eventually, these hypotheses lead to a system of 5 
differential equations: 



4 The 2D formulation 

A 2D treatment is essential for a reliable and physi- 
cal investigation of the effects of magnetic fields and 
turbulence on the solar interior, because the ID study 
does not allow for a full exploration of realistic mag- 
netic field configurations, and the 2D effects may be of 
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Fig. 1 The structure of a rotating 2D solar model is shown in terms of en, i.e. the variation of the radius compared to 
a non rotating ID, radially symmetric structure (left), the value of the colatitudinal acceleration gg (center), expressed in 
units of the solar surface gravity go, and the horizontal flux Fe (right). 
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The colatitudinal acceleration, Q 7 can be found via 
the horizontal component of the momentum equation 
and considering the Poisson equation: 



He + % + Tie tttt 

r oO 



(6) 



The 5 differential equations given above must be in- 
tegrated from the center to the surface, and from the 
pole to the equator. Provided that the structure is di- 
vided into N mass shell and M angles, the linearization 
procedure leads to 5(M-1)(N-1) equations for the 5MN 
unknown, i.e. the values of the 5 variables in the whole 
star. The system is closed with the following boundary 
conditions: 

• 3(M-1) central conditions, setting R=L=U=0 at the 
centre. 



• 2(M-1) boundary conditions at the surface, that ex- 
press the matching between the integration from the 
interior and the atmospheric treatment. 

• 5N conditions at the pole, i.e. that the 5 variables 
have null angular gradient at 9 = 0. 

The need of a high numerical accuracy demands one 
to find the shape of the equipotential surfaces before 
the integration of the above systems is performed. This 
is accomplished by using the condition of hydrostatic 
equilibrium, based on which we assume that pressure 
and its derivative with respect to mass are constant 
along the j > = co nst surfaces. The interested reader 



may find in lLi et al.l (|2006) and lLi et al.l (|2009f ) an ac- 
curate description of the numerical structure of the code 
used, and the techniques adopted to solve the full sys- 
tem of linearized equations. 



5 Applications of the 2D scheme 

The new formulation for the 2D approach was tested 
for two simple cases, namely i) a solar model where 
only the effects of a differential rotation in agreement 
with the observations of the solar spots is adopted, and 
ii) a structure where a purely toroidal magnetic field is 
allowed to vary in phase with the solar cycle. This step, 
considering the complexity of the numerical treatment 
requested, is mandatory before more complex magnetic 
field configurations, in which both the toroidal and the 
poloidal components are present, are investigated. 

5.1 The rotating model 

We studied the effects of differential rotation on the 
solar structure. No turbulence and magnetic field were 
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included. This requires one to neglect the 7i and T 
terms in eqs.1-5, and to consider only the 1Z vector, 
associated with rotation. 

The three panels of FigfT] show the modifications in- 
troduced by rotation in comparison with the standard 
ID radially symmetric solar model. We focus our at- 
tention on those quantities more sensitive to the 2D 
modeling, i.e. radius (left), colatitudinal acceleration Q 
(center), and transverse radiative flux, Fg (right). 

The radius variation is easily understood on the ba- 
sis of the expression for the centrifugal force, which in- 
creases from the pole to the equator, thus determining 
an expansion of the star at high colatitudes, and a con- 
traction at the poles; this is effectively what is observed 
in the left panel of Fig[TJ 

The contours for Q 7 shown in the central panel of 
FigE can be interpreted on the basis of eq.6 giving the 
expression for Q. Given the absence of the magnetic 
and turbulent contributions, the value of Q is given 
by the balance between the centrifugal horizontal term, 
\Q 2 r sin 20, and the derivative of pressure with respect 
to 9 calculated along the constant radius surfaces. This 
latter quantity is always positive, because the expan- 
sion of the structure towards the equator makes the 
r = const surfaces correspond to higher $ surfaces 
(hence, higher pressure) at higher 6. The two terms 
are therefore of opposite sign, with the centrifugal term 
dominating close to the equator, where the centrifugal 
acceleration is higher. Q is thus expected to be positive 



close to the equator. The pressure term is the dominant 
term close to the pole, where, indeed, we see from the 
central panel of Fig[T] that Q is negative. Note that Q 
vanishes as the equator is approached, since we expect 
purely radial contributions of all quantities there. 

The profiles of Fg are associated with the pressure 
and temperature profiles along the surface at constant 
radius. Like the pressure, even the temperature is ex- 
pected to increase moving at constant radius from the 
pole to the equator; this determines a flux in the op- 
posite direction, i.e. from high- to \ow-6 regions, which 
corresponds to a negative horizontal flux, as confirmed 
by the right panel of Fig[TJ Again, we may comment 
that the effects tend to vanish at the pole and at the 
equator and in the convective mantle, which is less sen- 
sitive to temperature variation than the radiative inte- 
rior. 

The luminosity (not shown) changes accordingly to 
what is discussed above. We stress the increase/decrease| 
of the radial flux on the equator/pole, which can be di- 
rectly associated to the behavior of the radius. The 
nuclear core, as expected, is almost insensitive to rota- 
tion. 

5.2 Toroidal field 

We started to investigate the effects of a magnetic field 
with a simple configuration, namely a toroidal field with 
two torus tubes that are parallel to the equatorial plane. 



g T (10~ 5 ) 



F fl /F (IP" 2 ) 



-0.97 -0.62 -0.28 0.07 0.42 0.77 



.26 -0.90 -0.53 -0.16 0.21 0.57 0.94 




0.0 



OA 



1 .0 



0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 

R/R R/R 

Fig. 3 Left: the same as Fig[3] but indicating the relative variation of temperature; right: the horizontal component of 
the radiative flux, Fg. 



The field B is assumed to peak in a region close to 
the bottom of the convective envelope, at a latitude 
of 9 = 30°, and to decay gradually with mass, with a 
gaussian profile. The region where B is located can be 
seen in the left panel of FigJ^l showing the ratio between 
the total pressure (including the magnetic contribution) 
and the pressure of the standard ID model, in which 
rotation, magnetic field and turbulence are neglected. 

The right panels of Fig[2] shows the relative varia- 
tion of density. The density contours confirm that the 
region where the magnetic perturbation is present ex- 
pands, due to the hydrostatic equilibrium condition, 
which demands a drop of the gas pressure. Note that, 
for continuity reasons, and because of the non-linearity 
of the stellar structure, wings of these effects reach also 
zones far away from the perturbation (i.e. the polar 
regions in the envelope). 

The temperature (shown in the left panel of Fig|3j) , 
is consequently decreased in the low-density regions, 
although the surface zone immediately above the per- 
turbation is indeed hotter, as demanded by hydrostatic 
equilibrium and the drop of density over there. The 
region where the perturbation is located is thus char- 
acterized by a strong horizontal temperature gradient, 
which is positive (i.e. temperature increases with the 
colatitude 9) in the higher 9 regions, and becomes neg- 
ative from the peak of the perturbation towards the 
equator. This consideration allows for a straightfor- 
ward understanding of the right panel of FigJ21 which 



shows a positive transverse flux in the upper part of 
the region encompassed by the magnetic field, and a 
negative Fg at lower latitudes. 



6 Conclusions 

A full comprehension of the mechanism driving the 
long-term variation of the total solar irradiance is far 
from being reached. The debate regarding the possible 
role played by periodic changes of the internal structure 
is still open, as confirmed by the many contributions in 
favor of a completely external solution, and other works 
stressing the importance to structural changes, most 
probably associated to a variable large-scale magnetic 
field. 

To understand how the presence of a magnetic per- 
turbation could contribute to the cyclic variation of the 
TSI, we developed a 2D solar evolution code, the two in- 
dependent coordinates being the mass enclosed within 
a given equipotential surface, m, and the colatitude, 
9. The full formulation consists in 5 differential equa- 
tions, to be integrated from the center to the surface, 
and from the pole to the equator. 

We applied this scheme to a test of the effects of a 
rotating solar model (note that both uniform and differ- 
ential rotation can be straightforwardly implemented), 
and to a structure accounting for the presence of a 
toroidal magnetic field. 
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More general configurations of the magnetic field will 
be investigated in future work. 
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